%%%%%%%%%%%
% fit_plr.m
%%%%%%%%%%%

% Data comes in as 2 x N matrix
% (i.e - rows are dimensions, columns are examples)
% so transpose to get as N x 2
Data = Data';

% Strip out any invalid (-1) rows
valid_idx = find(Data(:,1) ~= -1);
Data = Data(valid_idx, :);

[N,d] = size(Data);

x = Data(:,1);
y = Data(:,2);

Mean_x = mean(x);
Mean_y = mean(y);

Sxx = sum((x-repmat(Mean_x,N,1)).^2);
Sxy = sum((x-repmat(Mean_x,N,1)).*(y-repmat(Mean_y,N,1)));

% gradient
Grad = Sxy/Sxx;

% intercept
Offset = Mean_y - (Grad * Mean_x);

% predictions for data points
Pred_y = Grad*x + repmat(Offset,N,1);

% corresponding noise variance/MSE
Noise_var = sum((y-Pred_y).^2)/(N-2);

% 99th percentile of t-distribution with N-2 dof 
T_stat = tq(0.99,N-2);

% Save data
if mean(x) > 300
  save QQ_confs;
else
  save QT_confs;
end
